The Journal of Chemical Physics
● AIP Publishing
Preprints posted in the last 90 days, ranked by how well they match The Journal of Chemical Physics's content profile, based on 56 papers previously published here. The average preprint has a 0.02% match score for this journal, so anything above that is already an above-average fit.
Miao, Y.; Sastry, S.; Wang, J.; Mehmood, A.; Kim, M. T.-j.
Show abstract
It is valuable to calculate alchemical free energy changes in drug discovery and development. Thermodynamics Integration (TI) has been widely used in computational chemistry for estimating free energy changes with alchemical transformations. However, TI based on usually short Molecular Dynamics (MD) simulations often suffers from insufficient conformational sampling. Here, we have integrated Gaussian accelerated MD and TI (GaMD-TI) to enhance the conformational sampling and improve accuracy of free energy calculations. GaMD-TI has been demonstrated in model systems of alchemical changes in the Valine dipeptide and mutation cycle of the Alanine <-> Valine <-> Isoleucine (AVI) residues. Simulations showed that when GaMD boost potentials followed near-Gaussian distribution, the free energy change could be reweighted accurately through generalized cumulant expansion to the second order. The total free energy change often exhibited faster convergence using Selective GaMD (SGaMD) than using conventional MD (cMD). Accuracy of the free energy estimates from SGaMD-TI simulations was similar to or higher than those from cMD-TI simulations, although the differences were subtle for these small model systems. Meanwhile, dihedral angles in the model systems underwent significantly more frequent conformational transitions in SGaMD than in cMD, indicating improved sampling. Future studies are planned on larger systems with more complicated alchemical changes, such as ligand binding to proteins/nucleic acids and mutations at biomolecular binding interfaces. GaMD-TI should be broadly applicable to alchemical free energy calculations and therapeutic design.
Krott, L. B.; Puccinelli, T.; Oliveira, W. d.; Gomes, M. E. N.; Lomba, E.; Piazza, F.; Bordin, J. R.
Show abstract
Synapsin-1 is a multidomain neuronal protein containing extensive intrinsically disordered regions and is a key component of synaptic-vesicle condensates. Direct residue-level simulation of the collective organization of thousands of synapsin molecules remains computationally demanding. Here, we develop a coarse-grained description that connects residue-level CALVADOS 3 simulations to a one-particle-per-protein model. A potential of mean force between two synapsin molecules is obtained by umbrella sampling and represented by an isotropic effective interaction containing a short-range attractive region and a weak outer repulsive contribution. We compare two treatments of this interaction that differ only in the retention of the outer tail. Langevin dynamics simulations of effective proteins show aggregation upon cooling and compression in both models, but with markedly different collective organization. The shorter-ranged model progressively coarsens toward a single dense domain, whereas retaining the outer repulsive contribution favors the persistence of multiple mesoscale aggregates. The two models also display distinct relationships between aggregate size and particle mobility at low temperature. These results show that weak features of an effective protein-protein interaction can have pronounced consequences for collective synapsin organization at mesoscopic scales. TOC Graphic O_FIG O_LINKSMALLFIG WIDTH=200 HEIGHT=95 SRC="FIGDIR/small/745062v1_ufig1.gif" ALT="Figure 1"> View larger version (23K): org.highwire.dtl.DTLVardef@854288org.highwire.dtl.DTLVardef@d2fc52org.highwire.dtl.DTLVardef@1b37d5aorg.highwire.dtl.DTLVardef@eac583_HPS_FORMAT_FIGEXP M_FIG C_FIG
Gopi, S.; Qin, H.; Best, R. B.; Schuler, B.
Show abstract
Residue-level coarse-grained simulations provide a powerful route for modeling biomolecular condensates over length and time scales that are difficult to access with atomistic molecular dynamics. Coarse-grained models have been shown to reproduce many aspects of equilibrium phase behavior. However, it remains unclear to what extent such models can reproduce the relative timescales of molecular dynamics. Here, we examine this question for complex coacervates with markedly different dynamics, formed by the highly acidic intrinsically disordered protein prothymosin with four cationic partners: linker histone H1, protamine, polylysine, and polyarginine. Coexistence simulations using a residue-level coarse-grained model reproduce key equilibrium observables from experiments, including dense-phase concentrations, ionic-strength-dependent phase behavior, and chain dimensions in the dense and dilute phases. Dynamics are accelerated in these simulations, but a composition-specific time-rescaling factor captures the ionic-strength dependence of chain reconfiguration times within a given complex coacervate. In contrast, time rescaling is not transferable between dense and dilute phases or across condensate compositions and can depend on the chosen observable. These results show that agreement with measured equilibrium observables does not imply a universally transferable timescale for conformational dynamics in residue-level coarse-grained simulations. However, we find that the required time rescaling strongly correlates with the interaction energy of the protein chains, suggesting that the missing frictional effects arise from protein-protein interactions rather than solely from protein-solvent interactions, reminiscent of internal friction. Our findings highlight the need to combine thermodynamic validation with kinetic calibration when interpreting chain relaxation, molecular diffusion, and material properties from residue-level coarse-grained simulations of biomolecular condensates.
Nadeem, H.; Kleiman, D. E.; Shukla, D.
Show abstract
Adaptive sampling accelerates the exploration of conformational space in molecular dynamics (MD) simulations by repeatedly analyzing the accumulated trajectories and seeding a new round of simulations from informative configurations. A growing collection of adaptive sampling policies has been proposed, each built around a particular notion of what makes a configuration informative, yet these methods are scattered across separate and often incompatible implementations, which complicates their systematic comparison and their combined use in meta adaptive sampling schemes. Here, we present AdaptivePy, a compact and extensible Python framework that implements nine seed-selection policies behind a single configuration-driven interface, spanning simple population-based baselines, several established machine-learning and geometry-based methods, and two ensemble or meta sampling policies introduced in this work. We show that the shared implementation reproduces the characteristic selection behavior of each policy on a series of analytic benchmark landscapes. We also introduce a new adaptive sampling scheme that employs TS-DAR, a deep learning framework originally designed to identify transition states, into an acquisition criterion that drives the discovery of an entire multi-basin landscape starting from a single basin. We further demonstrate that the common interface enables meta adaptive sampling policies, which aggregate the rankings of several policies into a single set of seeds. AdaptivePy thereby provides a unified testbed for the adoption, benchmarking, and continued development of adaptive sampling methods for biomolecular MD simulations.
Dhibar, S.; Jana, B.
Show abstract
The process of drug unbinding is of immense importance in the field of biophysics and therapeutics. The behavior of these systems is greatly influenced by their thermodynamic and kinetic properties. Therefore, it is crucial to accurately estimate the ligand binding free energies and rate of ligand dissociation, yet these processes are often governed by rare event transitions that lie beyond the reach of standard brute-force molecular dynamics simulations. While enhanced sampling simulations offer a solution, their efficacy is strictly contingent upon the selection of appropriate collective variables (CVs) which is non-trivial for complex systems like protein-ligand complexes. In this study, we present a method to derive optimized CV from transition state region (TS) via an interpretable machine learning (ML) model, Elastic Net. By employing some physically intuitive order parameters, the derived optimized CV from the TS-region greatly accelerate ligand binding-unbinding transitions and achieves rapid free energy surface (FES) convergence across diverse systems including buried and solvent exposed active sites such as Trpsin-benzamidine complex, host-guest systems and sodium epoxidase etc. Intriguingly significant contribution of the ligand hydration is found in the optimized CV which depicts crucial role of solvent in driving ligand binding-unbinding transitions. The estimated binding free energies for different protein-ligand complexes match quite well with experiments, while maintaining a low computational cost. The derived optimized CV is also used to calculate the ligand residence times across different systems and calculated residence times are within the experimental range for all systems, again with very little computational costs. Moreover, we show that the optimized CV constructed from TS region via an interpretable ML model is transferable across diverse systems, offering a robust and scalable framework for drug discovery and investigation of complex biomolecular recognition.
Akram, S. A.; Liu, J.; Chattaraj, A.; Gautam, S.; Chakrabarti, B.; Shakhnovich, E.; Schmit, J. D.
Show abstract
Phase separating biopolymers are often described by the sticker-spacer model, which accounts for non-attractive linkers that connect discrete attractive moieties. "Associative polymer" theories describe these systems as two coexisting fluids of stickers and spacers, respectively. In such theories the stickers and spacers are decoupled in the sense that sticker binding does not affect spacer conformations or vice versa. Informed by coarse-grained simulations, we develop a phenomenological mean-field theory to determine the microscopic coupling between the sticker and spacer degrees of freedom. Accounting for linkers connecting adjacent spacers introduces free energy contributions due to conformational entropy of the spacer and the binding combinatorics of the stickers. These contributions have opposite dependencies on spacer length because longer spacers allow a greater volume to search for binding partners but incur a greater entropy cost upon bond formation. The length dependent terms from these contributions cancel leaving a length-independent elastic term arising from the fact that a small amount of spacer stretching lowers the free energy by increasing the number of possible binding partners. The cancellation of spacer length from the theory justifies the decoupling approximation and supports mean-field treatments of heterogeneous spacers.
Elgendy, A.; Zeipelt, A. P.; Schäfer, L. V.
Show abstract
Molecular dynamics (MD) simulations of slow biomolecular processes, such as exploration of the conformational ensembles of intrinsically disordered proteins (IDPs), are computationally demanding. Although coarse-grained (CG) models can substantially speed up the simulations compared to all-atom MD, the sampling challenge can still be significant for large systems and long time scales. Here, we present light Martini water, a low-viscosity water model that accelerates sampling in MD simulations with the Martini CG force field. We systematically reduced the mass of the Martini water beads and verified stable, accurate integration of the equations of motion with 20 fs time steps, as typically used in Martini simulations. Light Martini water has a reduced mass of 20 amu (compared to 72 amu in the standard water model), yielding up to a 2.68-fold increase in the sampling rate of IDP chain reconfiguration in water and a 16 % increase in the lateral diffusion of lipids in a POPC bilayer. Equilibrium properties remained unaffected by the mass scaling, and the speedup was achieved without compromising simulation accuracy. The water model is trivial to implement, has no computational overhead, and should be universally applicable to Martini simulations.
Yang, Y.
Show abstract
Gaussian accelerated molecular dynamics (GaMD) enhances conformational sampling by adding a smooth boost potential without requiring predefined collective variables, but an engine-integrated implementation has not been available in GROMACS. Here, we implement total-, dihedral-, and dual-boost GaMD in GROMACS 2025.4, including staged energy-statistics collection, GPU-based bias evaluation and force scaling, restart support, and outputs required for cumulant-based free-energy reweighting. The implementation was evaluated using four benchmark systems spanning conformational free energies, protein folding, and ligand recognition. For alanine dipeptide, a reweighted 100 ns GaMD trajectory recovered the major free-energy basins and rotational barriers in overall agreement with a 1000 ns conventional MD simulation. For chignolin and TC5b, all three independent trajectories for each system sampled native-like folded states from extended conformations within 300 ns and 1 s, respectively; the best TC5b structure had a minimum backbone RMSD of 0.03 nm from the experimental structure. In the benzene-T4 lysozyme system, two of five independent 500 ns trajectories captured both ligand binding and dissociation, yielding a bound pose with a minimum ligand RMSD of 0.06 nm from the crystal structure. Across all four systems, the boost-potential distributions were approximately Gaussian, and second-order cumulant reweighting resolved the expected conformational and binding free-energy basins. These results demonstrate that GROMACS-GaMD provides a practical, GPU-enabled, collective-variable-free enhanced-sampling framework for biomolecular free-energy calculations, protein folding, and ligand-binding studies.
Ghermezcheshme, H.; Bisht, K.; Yesylevskyy, S. O.; Khandelia, H.
Show abstract
Nanoplastic interactions with cell membranes can influence particle uptake, translocation, and biological effects, making their molecular description important for assessing the environmental and biological consequences of plastic pollution. However, molecular simulations generally represent nanoplastics using monodisperse polymer chains, whereas nanoplastics formed by environmental degradation are polydisperse. Here, we construct 4 nm polyethylene terephthalate (PET), polyethylene (PE), and polystyrene (PS) nanoparticles with lognormal molecular-weight distributions reported for degraded polymers, and with melt-like chain packing and entanglement. Using Martini coarse-grained molecular dynamics, we simulate their interactions for 10 microseconds with three lipid bilayers of increasing complexity. PE and PS insert rapidly and remain largely intact, whereas PET remains surface-associated and progressively releases individual chains. Fluid, polyunsaturated lipid environments undergo greater deformation and promote substantially more chain release, up to 56% of the chains for PET and 43% for PS. Release begins with the shortest chains and progressively extends to longer chains, while semicrystalline PE remains largely intact. In contrast, a size-matched monodisperse PS particle releases no chains and minimally perturbs the membrane. These results identify the molecular-weight distribution as a key determinant of nanoplastic-membrane interactions, and indicate that uniform model particles underestimate the release of polymer chains into biological membranes.
Argun, B. R.; Stachowiak, J.; Ren, P.
Show abstract
Recent experiments show that protein condensates sitting on opposite surfaces of a flat lipid membrane move together and prefer to overlap, even though they cannot touch each other. This points to an indirect, membrane-mediated interaction. Two mechanisms could be responsible: a curvature-induced interaction, which is energetic in origin, and a fluctuation-induced interaction, which is entropic. Here we study both with coarse-grained molecular dynamics simulations, using Cookes implicit-solvent lipid model together with a generic bead-spring polymer model for the condensate. We compute the potential of mean force between two condensates across the membrane. For condensates of the same size, full overlap is unfavorable, and the pair instead settles into a partially overlapping state that bends the membrane into an S-like shape. When the two condensates differ strongly in size, full overlap becomes favorable. We explain this with a simple geometric picture. The condensate wets the membrane as a thin film and imposes curvature only along its rim, while membrane tension flattens the membrane under its interior. The resulting ring of curvature can trap a smaller condensate on the opposite side. We also compare the bending undulations and the effective bending modulus of a bare membrane, a membrane with one condensate, and a membrane with condensates on both sides. A wetting condensate suppresses the undulation modes and stiffens the membrane, but whether this makes overlap entropically favorable remains inconclusive. Our results indicate that the coupling is driven mainly by curvature, and that it depends on the wetting mechanism and on the membrane tension.
ADUPA, V.; Polet, J. D.; Dekker, M.; Onck, P. R.
Show abstract
Polyglutamine (polyQ) aggregation plays a central role in several neurodegenerative diseases, including Huntington's disease. DNAJB6b, a molecular chaperone involved in protein quality control, is known to efficiently suppress polyQ aggregation, but its anti-aggregation mechanism remains unclear. In this work we investigate the interaction between DNAJB6b and the polyQ region (Q48) of mutant Huntingtin Exon 1 (mHttEx1) using a custom-built coarse-grained molecular dynamics model. The model incorporates a two-bead-per-amino-acid representation with hydrogen bonding (termed 2BPA-HB), and is calibrated against all-atom molecular dynamics data in terms of geometry, hydrophobicity, and hydrogen bonding. The model reproduces the tertiary structure of DNAJB6b and its interactions with Q48, and reveals an inverse correlation between DNAJB6b concentration and Q48 aggregation propensity. Our simulations show that DNAJB6b co-condensates with polyQ molecules, thereby shielding the polyQ from forming the intermolecular hydrogen bonds necessary for amyloid formation. The 2BPA-HB CGMD model en- ables efficient exploration of DNAJB6b conformations, supporting future studies of chaperone-mediated aggregation suppression and therapeutic development.
Paesani, M.; Ilie, I. M.
Show abstract
Successful intracellular delivery requires a nanocarrier to adhere to the cell membrane, be wrapped and internalized and finally released into the cytoplasm before it is degraded. While the initial steps have been intensively studied, the intracellular release of a nanocarrier has received comparatively less attention and remains poorly understood. To understand the physical mechanisms that govern intracellular release, we investigate the escape of nanocarriers (radii 4.5{sigma} and 7{sigma} ) from a lipid vesicle using coarse-grained simulations. We represent a nanocarrier using the metaparticle model and the membrane using the three-beads Cooke-Deserno model. Starting from a fully endocytosed state, we systematically vary the repulsion between the nanoparticle and the lipid head groups and compare three release scenarios: passive release, localized inside-out activity from a drug-like source confined within the particle and activity distributed over the entire nanoparticle. An elastic-energy analysis shows that in a tensionless membrane escape requires crossing a large pore-expansion barrier of {approx} 68kBT for the smaller particle and {approx} 110kBT for the larger, making passive release thermodynamically improbable. Our results show that the smaller particle nonetheless escapes above a threshold repulsion, whereas the larger one remains trapped in a metastable state in which local membrane pores form but fail to grow beyond the particle diameter to allow successful release. Adding activity lowers the nucleation barrier and promotes earlier pore opening, but its effectiveness depends on how it is applied, i.e., a localized internal source assists release of the smaller particle yet is insufficient to systematically free the larger one, while activity distributed across the whole nanoparticle drives complete release even in the larger case. These results show that productive escape is governed by whether the resulting pore can evolve into a geometrically accessible release pathway, providing physical design principles for intracellular delivery by nanocarriers
Premathilaka, M. B.; Hansmann, U. H. E.
Show abstract
Motivated by correlations between SARS-COV-2 infections and transthyretin (TTR) amyloidosis we use molecular dynamics simulations to study whether three protein fragments from the SARS-CoV-2 virus can initiate dissociation of the transthyretin (TTR) homotetramer as the first step on the pathway to TTR amyloidosis. We find that the peptides reduce the frequency of transient hydrogen bonds at the monomer-monomer interfaces, decreasing in this way the stability of the TTR tetramer.
Zhou, M.; Jia, C.
Show abstract
Moment-closure methods provide a widely used approach for approximating low-order statistics of stochastic biochemical reaction networks, but the resulting closed equations may exhibit physically meaningless dynamics. Here, we systematically investigate the stability of second-order moment-closure approximations for a bursty gene expression model with nonlinear protein degradation. We consider the Gaussian, log-normal, and negative-binomial closures. By analyzing the resulting dynamical systems, we show that the Gaussian closure fails the proposed stability criterion for all parameter values, whereas both the log-normal and negative-binomial closures possess a unique physically accessible fixed point and are stable if and only if kB(1 + 2B) [≥] 1, where B is the mean burst size and k is the dimensionless protein synthesis rate. Moreover, this condition is automatically satisfied when the steady-state mean protein copy number is at least one, suggesting that it is mild in biologically relevant regimes. Numerical comparisons over a broad parameter range further show that the negative-binomial closure is the most accurate, while the Gaussian closure is the least accurate. We further develop higher-order negative-binomial moment-closure schemes and demonstrate that increasing the closure order substantially improves predictive accuracy. These results highlight the importance of considering dynamical stability, in addition to numerical accuracy, when assessing moment-closure approximations.
Tulli, L.; Lodola, A.; Mulholland, A. J.; Oliveira, A. S. F.
Show abstract
The epidermal growth factor receptor (EGFR) is a tyrosine kinase receptor that plays a fundamental role in regulating cellular proliferation, survival, and differentiation. EGFR's kinase domain catalyses the autophosphorylation that drives downstream signalling and is a known hotspot for cancer-associated mutations. Such substitutions often enhance receptor activation, driving oncogenesis and therapeutic resistance, with major implications for prognosis and treatment choice; consequently, a detailed understanding of this domain's functional dynamics, and how mutations reshape them, is essential. Here, we apply dynamical-nonequilibrium molecular dynamics (D-NEMD) simulations to probe how changes at the ATP-binding site are transmitted throughout the kinase domain, ultimately shaping the behaviour of distal regions in the wild-type protein and in two of the most clinically observed non-small cell lung cancer variants, L858R and L858R+T790M. D-NEMD reveals previously uncharacterised communication networks linking the ATP pocket to key functional motifs, including the P-loop, C-helix, hinge region, activation loop, and F-helix, showing the intricate allosteric connectivity within this domain. Notably, numerous cancer-related mutation sites, despite being distant from the ATP site, map onto or lie close to these communication pathways, indicating that D-NEMD simulations can identify functionally relevant allosteric positions and potential mutation sites in EGFR. Our simulations further reveal that the two variants rewire EGFR's internal signal propagation network in distinct ways: L858R diminishes the transmission of structural changes to distal regions of the C-lobe, whereas introducing T790M in addition to L858R partially restores wild-type-like behaviour through compensatory rerouting of early-stage dynamical pathways.
Zhang, G.; Leng, Y.
Show abstract
Interpreting single-molecule force spectroscopy (SMFS) data from optical tweezers with steered molecular dynamics (SMD) is hindered by a severe timescale mismatch: SMD trajectories span nanoseconds to microseconds, while SMFS experiments unfold over milliseconds to seconds, forcing SMD to pull six to nine orders of magnitude faster than experiment and systematically distorting rupture forces, molecular pathways, and kinetic intermediates. We introduce force-probe coupled steered molecular dynamics (fp-SMD), which explicitly propagates the force probe (the optically trapped bead) alongside the biomolecule within a single classical Hamiltonian description of the coupled molecule-linker-bead system. A Predict-Correct Trajectory Propagation (PCTP) scheme resolves the resulting computational bottleneck by decoupling the fast molecular relaxation from the slow, macroscopic bead motion, using metadynamics and transition-state theory to identify the kinetically resolved molecular intermediate reached at each macroscopic time step. On the cholesterol-{beta}-cyclodextrin ({beta}-CD) host-guest complex, fp-SMD reproduces the experimental rupture-force histogram directly at matched conditions, without extrapolation from standard SMD. fp-SMD further resolves how applied force reshapes the molecular escape pathway itself--cholesterol dissociates via a high-angle route at low force and a direct, small-angle route at high force--providing, to our knowledge, the first direct structural evidence for force-induced pathway switching within small host-guest systems believed to underlie catch-slip bond behavior. fp-SMD/PCTP thus offers a general route to simulating force-probe-coupled biomolecular mechanics at experimentally relevant timescales.
Zhang, Y.; Sood, A.; Athreya, A.; Zhang, B.
Show abstract
Biomolecular condensates formed by intrinsically disordered proteins (IDPs) are often described using stickers-and-spacers models, in which specific sticker motifs form reversible crosslinks and spacer regions modulate phase behavior. Recent theory predicts that heterogeneous nonspecific spacer interactions can promote condensation but may also disrupt sticker-mediated organization. Here, we develop an off-lattice coarse-grained stickers-and-spacers polymer model for continuous molecular dynamics simulations and implement it in the GPU-accelerated OpenABC package. The model uses a directional sticker-sticker interaction to encode limited valency through interaction geometry, producing effectively one-to-one sticker binding without explicit bond assignment. Simulations of one-component systems show that sticker affinity and multivalency promote porous, network-like condensates, while nonspecific spacer interactions can also drive phase separation but produce more compact, spacer-dominated dense phases. When both interaction types are present, strong spacer heterogeneity reduces sticker conversion, suppresses sticker mobility, and disrupts the sticker-mediated network. In two-component systems, specific sticker interactions buffer client recruitment into host condensates, while nonspecific spacer interactions produce reservoir-dependent uptake. These results support a tradeoff in which spacer heterogeneity promotes condensation at the cost of condensate organization and compositional robustness, providing a physical rationale for the suppression of promiscuous spacer interactions in low-complexity IDP regions.
Chan, B.; Rubinstein, M.
Show abstract
In the active loop extrusion model, the cohesin protein complex creates chromatin loops in eukaryotic cells. Extrusion maintains topologically associated domains (TADs), which are contiguous segments of chromatin that preferentially colocalize in space and are typically bounded by CTCF proteins that pause cohesin translocation. Here, we model active loop extrusion with hybrid molecular dynamics - Monte Carlo simulations in entangled flexible linear polymer melts. Intra-chain contact probabilities of polymers with active loop extrusion are enhanced compared to their equilibrium, passive counterparts. Extrusion causes the size of chain segments to be much smaller than in passive melts. While the overlap parameter in passive melts without extrusion monotonically increases with segment length, it is nonmonotonic in active melts and on the order of unity within the parameters of this study. Active loop extrusion suppresses contacts between TADs in favor of intra-TAD contacts. Reduction of overlaps between chain segments dilutes entanglements in active melts. Depending on parameters, active extrusion without TADs may induce more compact conformations than with TADs, due in part to fractal loopy globule-like dynamics. This work suggests that active loop extrusion reduces overlaps between TADs, contributing to effective gene regulation by cis-regulatory elements.
Herb, N.; Brajkovic, M.; DArrigo, G.; Kokh, D. B.; Wade, R. C.
Show abstract
Interleukin-13 (IL-13) is an immunomodulatory cell signaling cytokine that has been implicated in neurodegenerative disease and chronic inflammation. IL-13 binds to its low and high affinity receptors, IL-13 receptor 1 (IL-13R1) and IL-13 receptor 2 (IL-13R2), respectively, with residence times that vary accordingly. As the binding kinetics of the cytokine-receptor complexes influence cellular responses, we employed the molecular dynamics (MD) simulation-based{tau} -random acceleration molecular dynamics method ({tau}RAMD) to compute relative residence times for wild-type (WT) IL-13 and 19 IL-13 mutants to the two receptors. Comparison with experimental kinetic data shows that the{tau} RAMD computations capture the trends in residence times. Analysis of simulated dissociation trajectories of the cytokine-receptor complexes reveals two distinct dissociation pathways of IL-13 from each of the receptors. This study thus pinpoints key determinants of the interaction of IL-13 with its receptors which could be targeted for therapeutic design. Statement of SignificanceCytokines are regulatory proteins that bind to cell surface receptors and thereby send signals to the cellular interior. Interleukin-13 (IL-13) is a cytokine that has a low and a high affinity receptor. It has important physiological roles, and its deregulation is involved in diseases such as atopic dermatitis and asthma. We employed a molecular dynamics simulation-based method to compute the effects of changes in the sequence of IL-13 on the lifetimes of complexes of IL-13 and its receptors. Comparison with experiments supports the validity of the computational approach and analysis of the simulations reveals two distinct ways in which IL-13 dissociates from each receptor. These results thus provide a map for targeting IL-13 - receptor interactions for the design of therapeutics.
Bhandari, T. R.; Kremer, K.; Girard, M.
Show abstract
Simulations have emerged as a pillar in biophysics to understand behavior at the molecular scale, in particular for proteins such as intrinsically disordered proteins. These often show transient folding with long-lived states that are challenging to efficiently sample using conventional simulations. Here, we show that Hamiltonian replica-exchange at the level of hydrogen bonds can be used to systematically explore the conformational landscape of proteins. Investigating the IM30 protein, a partially alpha helical protein, under the PLUM coarse-grain force-field, we show that this approach straightforwardly allows characterization of local helix stability, which would be otherwise hard to obtain. For instance, we show that some regions of the conserved coiled-coil structure are less stable than others. By combining the different replicas with principal component analysis, we are able to dissect the conformational landscape of the IM30 protein. Taken together, this indicates that enhanced sampling remains an important part of simulations, allowing insights that cannot be easily obtained in traditional simulations.